clear all;
matlabpool(12);

%Setting Options
optionsLoose = optimset('TolX', 1e-5, 'TolFun',1e-5, 'MaxFunEvals', 1e5, 'MaxIter', 1e5);
optionsTight = optimset('TolX', 1e-8, 'TolFun',1e-8, 'MaxFunEvals', 1e8, 'MaxIter', 1e8);

UtilFuncParam = 1;
ChoiceProbForm = 1;
NumBootReps = 250;
RandomizationSeed = 12345;

%Loading Data
DATA = csvread('CCES_Matlab_Senate.csv', 1, 0);

VoteDem = DATA(:,1);
VoteDem_1 = DATA(:,1);
VoterIdeal = DATA(:,2);
RepPartyIdeal = DATA(:,3);
DemPartyIdeal = DATA(:,4);
RepCandIdeal = DATA(:,5);
DemCandIdeal = DATA(:,6);
IncRep = DATA(:,7);
IncDem = DATA(:,8);
ExpChallRep = DATA(:,9);
ExpChallDem = DATA(:,10);
Year2008 = DATA(:,11);
Year2010 = DATA(:,12);
Year2012 = DATA(:,13);
Age = DATA(:,14);
State_FE=DATA(:,15:63);
State_Index = DATA(:,65); 

CandTypeDem = horzcat(IncDem, ExpChallDem);
CandTypeRep = horzcat(IncRep, ExpChallRep);
clear DATA;


Lambdas = vertcat(0.5, 0.5);
RHS = [Year2008, Year2010, Year2012];
Betas0 = vertcat(1, 1, zeros(length(RHS(1,:)), 1)); Theta_0 = vertcat(Lambdas, Betas0, 1);

[thetaTightProbit1G, fvalProbit1G] = fminunc(@(Theta) IndivVoterObjFun_Senate(Theta,VoteDem, VoterIdeal, CandTypeDem, CandTypeRep, DemCandIdeal, DemPartyIdeal, RepCandIdeal, RepPartyIdeal, RHS, UtilFuncParam, ChoiceProbForm), Theta_0, optionsTight);

[SEsBoot1, thetaBoot1, exitflagBoot1] = getVarCovBootstrapCluster_Senate(thetaTightProbit1G, VoteDem, VoterIdeal, CandTypeDem, CandTypeRep, DemCandIdeal, DemPartyIdeal, RepCandIdeal, RepPartyIdeal, RHS, UtilFuncParam, ChoiceProbForm, NumBootReps, optionsLoose, optionsTight, RandomizationSeed, State_Index);

%%
[IncAdv_Mat_Probit1] = getIncAdv_Senate(thetaTightProbit1G,VoteDem, VoterIdeal, CandTypeDem, CandTypeRep, DemCandIdeal, DemPartyIdeal, RepCandIdeal, RepPartyIdeal, RHS, UtilFuncParam, 1,State_Index);
csvwrite('thetaTightProbit1_Senate.csv', thetaTightProbit1G);
csvwrite('fvalProbit1_Senate.csv', fvalProbit1G);
csvwrite('IncAdv_Mat_Probit1_Senate.csv', IncAdv_Mat_Probit1);
csvwrite('thetaBoot_Probit1_Senate.csv', thetaBoot1);

for b = 1:NumBootReps
    [IncAdv_Mat_Probit1_Boot(:,b)] = getIncAdv_Senate(thetaBoot1(:,b),VoteDem, VoterIdeal, CandTypeDem, CandTypeRep, DemCandIdeal, DemPartyIdeal, RepCandIdeal, RepPartyIdeal, RHS, UtilFuncParam, 1,State_Index);
end

IncAdv_Mat_Probit1_SE = sqrt((1/(NumBootReps-1)).*sum((IncAdv_Mat_Probit1_Boot - repmat(mean(IncAdv_Mat_Probit1_Boot,2), 1, NumBootReps)).^2,2));

csvwrite('IncAdv_Mat_Probit1_Boot_Senate.csv', IncAdv_Mat_Probit1_Boot);
csvwrite('IncAdv_Mat_Probit1_SE_Senate.csv', IncAdv_Mat_Probit1_SE);


%%
RHS = [Year2008, Year2010, Year2012, State_FE];
Betas0 = vertcat(1, 1, zeros(length(RHS(1,:)), 1)); Theta_0 = vertcat(Lambdas, Betas0, 1);
[thetaTightProbit2G, fvalProbit2G] = fminunc(@(Theta) IndivVoterObjFun_Senate(Theta,VoteDem, VoterIdeal, CandTypeDem, CandTypeRep, DemCandIdeal, DemPartyIdeal, RepCandIdeal, RepPartyIdeal, RHS, UtilFuncParam, ChoiceProbForm), Theta_0, optionsTight);
[SEsBoot2 thetaBoot2 exitflagBoot2] = getVarCovBootstrapCluster_Senate(thetaTightProbit2G, VoteDem, VoterIdeal, CandTypeDem, CandTypeRep, DemCandIdeal, DemPartyIdeal, RepCandIdeal, RepPartyIdeal, RHS, UtilFuncParam, ChoiceProbForm, NumBootReps, optionsLoose, optionsTight, RandomizationSeed, State_Index);

[IncAdv_Mat_Probit2] = getIncAdv_Senate(thetaTightProbit2G,VoteDem, VoterIdeal, CandTypeDem, CandTypeRep, DemCandIdeal, DemPartyIdeal, RepCandIdeal, RepPartyIdeal, RHS, UtilFuncParam, 1,State_Index);
csvwrite('thetaTightProbit2_Senate.csv', thetaTightProbit2G);
csvwrite('fvalProbit2_Senate.csv', fvalProbit2G);
csvwrite('IncAdv_Mat_Probit2_Senate.csv', IncAdv_Mat_Probit2);
csvwrite('thetaBoot_Probit2_Senate.csv', thetaBoot2);


for b = 1:NumBootReps
    [IncAdv_Mat_Probit2_Boot(:,b)] = getIncAdv_Senate(thetaBoot2(:,b),VoteDem, VoterIdeal, CandTypeDem, CandTypeRep, DemCandIdeal, DemPartyIdeal, RepCandIdeal, RepPartyIdeal, RHS, UtilFuncParam, 1,State_Index);
end

IncAdv_Mat_Probit2_SE = sqrt((1/(NumBootReps-1)).*sum((IncAdv_Mat_Probit2_Boot - repmat(mean(IncAdv_Mat_Probit2_Boot,2), 1, NumBootReps)).^2,2));

csvwrite('IncAdv_Mat_Probit2_Boot_Senate.csv', IncAdv_Mat_Probit2_Boot);
csvwrite('IncAdv_Mat_Probit2_SE_Senate.csv', IncAdv_Mat_Probit2_SE);


%%

DATA = csvread('CCES_Matlab_IndivCovariates_Senate.csv', 1, 0);

VoteDem = DATA(:,1);
VoterIdeal = DATA(:,2);
RepPartyIdeal = DATA(:,3);
DemPartyIdeal = DATA(:,4);
RepCandIdeal = DATA(:,5);
DemCandIdeal = DATA(:,6);
IncRep = DATA(:,7);
IncDem = DATA(:,8);
ExpChallRep = DATA(:,9);
ExpChallDem = DATA(:,10);
Year2008 = DATA(:,11);
Year2010 = DATA(:,12);
Year2012 = DATA(:,13);
Age = DATA(:,14);
State_FE=DATA(:,15:63);
State_Index = DATA(:,65); 

HighSchool = DATA(:,66);
TwoYearCollege = DATA(:,67);
FourYearCollege = DATA(:,68);
PostBA = DATA(:,69);
Female = DATA(:,70);
Af_Am = DATA(:,71);
Latino = DATA(:,72);
Asian = DATA(:,73);
AttendRelServiceWeek = DATA(:,74);
Income = DATA(:,75);


CandTypeDem = horzcat(IncDem, ExpChallDem);
CandTypeRep = horzcat(IncRep, ExpChallRep);

%%
RHS = [Year2008, Year2010, Year2012, State_FE, HighSchool, TwoYearCollege, FourYearCollege, PostBA, Female, Af_Am, Latino, Asian, AttendRelServiceWeek, Income, Age];


Betas0 = vertcat(1, 1, zeros(length(RHS(1,:)), 1)); Theta_0 = vertcat(Lambdas, Betas0, 1);
[thetaTightProbit4G, fvalProbit4G] = fminunc(@(Theta) IndivVoterObjFun_Senate(Theta,VoteDem, VoterIdeal, CandTypeDem, CandTypeRep, DemCandIdeal, DemPartyIdeal, RepCandIdeal, RepPartyIdeal, RHS, UtilFuncParam, ChoiceProbForm), Theta_0, optionsTight);
[SEsBoot4, thetaBoot4, exitflagBoot4] = getVarCovBootstrapCluster_Senate(thetaTightProbit4G, VoteDem, VoterIdeal, CandTypeDem, CandTypeRep, DemCandIdeal, DemPartyIdeal, RepCandIdeal, RepPartyIdeal, RHS, UtilFuncParam, ChoiceProbForm, NumBootReps, optionsLoose, optionsTight, RandomizationSeed, State_Index);

[IncAdv_Mat_Probit4, IncAdvExpChallDem_Weighted_Holder4, IncAdvExpChallRep_Weighted_Holder4] = getIncAdv_Senate(thetaTightProbit4G,VoteDem, VoterIdeal, CandTypeDem, CandTypeRep, DemCandIdeal, DemPartyIdeal, RepCandIdeal, RepPartyIdeal, RHS, UtilFuncParam, 1,State_Index);


%%
myfiguresize =  [100 1 10 8];
set(gcf, 'PaperPosition', myfiguresize);
fh = figure(1);
get(gca)
set(gcf, 'color', 'white');

hold 'on'
[f_Dem,xi_Dem] = ksdensity(IncAdvExpChallDem_Weighted_Holder4);
[f_Rep,xi_Rep] = ksdensity(IncAdvExpChallRep_Weighted_Holder4);
plot(xi_Dem, f_Dem, 'b--')
plot(xi_Rep, f_Rep, 'r')
legend('Democratic Inc Adv', 'Republican Inc Adv')
legend boxoff
title('Density of State Average Incumbency Advantage (Senate Elections)')
hold 'off'

csvwrite('IncAdvExpChallDem_Weighted_Holder4_Senate.csv', IncAdvExpChallDem_Weighted_Holder4);
csvwrite('IncAdvExpChallRep_Weighted_Holder4_Senate.csv', IncAdvExpChallRep_Weighted_Holder4);

saveas(fh, 'IncAdv_Density_Senate', 'epsc'); 
[h, p, ks2stat] = kstest2(IncAdvExpChallDem_Weighted_Holder4, IncAdvExpChallRep_Weighted_Holder4)
csvwrite('pValKSTest_Senate.csv', p)
close

csvwrite('thetaTightProbit4_Senate.csv', thetaTightProbit4G);
csvwrite('fvalProbit4_Senate.csv', fvalProbit4G);
csvwrite('IncAdv_Mat_Probit4_Senate.csv', IncAdv_Mat_Probit4);
csvwrite('thetaBoot_Probit4_Senate.csv', thetaBoot4);

for b = 1:NumBootReps
    [IncAdv_Mat_Probit4_Boot(:,b)] = getIncAdv_Senate(thetaBoot4(:,b),VoteDem, VoterIdeal, CandTypeDem, CandTypeRep, DemCandIdeal, DemPartyIdeal, RepCandIdeal, RepPartyIdeal, RHS, UtilFuncParam, 1,State_Index);
end

IncAdv_Mat_Probit4_SE = sqrt((1/(NumBootReps-1)).*sum((IncAdv_Mat_Probit4_Boot - repmat(mean(IncAdv_Mat_Probit4_Boot,2), 1, NumBootReps)).^2,2));

csvwrite('IncAdv_Mat_Probit4_Boot_Senate.csv', IncAdv_Mat_Probit4_Boot);
csvwrite('IncAdv_Mat_Probit4_SE_Senate.csv', IncAdv_Mat_Probit4_SE);


Polarization_Factor_Vec= [0.5:0.1:2]';

for k = 1:length(Polarization_Factor_Vec)
    [IncDem_Weighted_Polarization(k,:), IncRep_Weighted_Polarization(k,:),  ExpChalltoIncDem_Weighted_Polarization(k,:), ExpChalltoIncRep_Weighted_Polarization(k,:)] = getPolarization_Senate(thetaTightProbit4G,VoteDem, VoterIdeal, CandTypeDem, CandTypeRep, DemCandIdeal, DemPartyIdeal, RepCandIdeal, RepPartyIdeal, RHS, UtilFuncParam, 1,State_Index, Polarization_Factor_Vec(k,:));
end

%%
close

myfiguresize =  [100 1 10 8];
set(gcf, 'PaperPosition', myfiguresize);
fh = figure(1);
get(gca)
set(gcf, 'color', 'white');

hold 'on'
plot(Polarization_Factor_Vec, IncDem_Weighted_Polarization, '--bs', 'MarkerSize', 8, 'Marker', '+')
plot(Polarization_Factor_Vec, IncRep_Weighted_Polarization, '--rs', 'MarkerSize', 8)
ylim([-0.04 0.14])
ylabel('Estimated Incumbency Advantage')
xlabel('Polarization Multiplier')
title('Simulated Incumbency Advantage as Function of Party Polarization (Senate Elections)')
legend('Democrats (Incumbents to Challengers)', 'Republicans (Incumbents to Challengers)')
legend('Location', 'NorthWest')
legend boxoff
hold 'off'


saveas(fh, 'Polarization_Senate', 'epsc'); 
csvwrite('IncDem_Weighted_Polarization_Senate.csv', IncDem_Weighted_Polarization);
csvwrite('IncRep_Weighted_Polarization_Senate.csv', IncRep_Weighted_Polarization);
close

LogLikelihoodNoImpute = horzcat(-fvalProbit1G, -fvalProbit2G, -fvalProbit4G);

csvwrite('LogLikelihoodNoImputeCollect_Probit_Senate.csv', LogLikelihoodNoImpute);

labels = vertcat(char('Party Weight (Incumbent)',  'Party Weight (Chall)', 'Incumbent', 'Challenger'));
blanks = char('');

FID = fopen('CCESTable1_Probit_Senate.tex', 'w');
fprintf(FID, '\\begin{table}[t]\\centering \n');
fprintf(FID,  '\\begin{tabular}{l*{3}{c}} \\hline\\hline \n');
fprintf(FID, '&\\multicolumn{1}{c}{(1)}& \\multicolumn{1}{c}{(2)}&\\multicolumn{1}{c}{(3)}&\\\\ \n');
fprintf(FID, '&\\multicolumn{1}{c}{}&\\multicolumn{1}{c}{}&\\multicolumn{1}{c}{}&\\\\ \n');
fprintf(FID, '\\hline \n');
for k=1:4
     fprintf(FID, '%s & %8.3f & %8.3f & %8.3f \\\\ ', labels(k,:), round(thetaTightProbit1G(k)*10000)/10000, round(thetaTightProbit2G(k)*10000)/10000, round(thetaTightProbit4G(k)*10000)/10000);
     fprintf(FID, '%s & %s & %s & %s \\\\ \n', blanks, strcat('(', num2str(SEsBoot1(k)), ')'), strcat('(', num2str(SEsBoot2(k)), ')'), strcat('(', num2str(SEsBoot4(k)), ')'));   
end
             fprintf(FID, '\\hline \n');
             
labels_IncAdv = vertcat(char('Dem Inc Adv',  'Rep Inc Adv', 'Dem Signaling (\%)', 'Rep Signaling (\%)'));
             
for k=1:2
     fprintf(FID, '%s & %8.3f & %8.3f & %8.3f \\\\ ', labels_IncAdv(k,:), round(IncAdv_Mat_Probit1(k)*10000)/10000, round(IncAdv_Mat_Probit2(k)*10000)/10000, round(IncAdv_Mat_Probit4(k)*10000)/10000);
     fprintf(FID, '%s & %s & %s & %s \\\\ \n', blanks, strcat('(', num2str(IncAdv_Mat_Probit1_SE(k)), ')'), strcat('(', num2str(IncAdv_Mat_Probit2_SE(k)), ')'), strcat('(', num2str(IncAdv_Mat_Probit4_SE(k)), ')'));   
end

for k=3:4
     fprintf(FID, '%s & %8.3f & %8.3f & %8.3f \\\\ ', labels_IncAdv(k,:), round(IncAdv_Mat_Probit1(k+2)*10000)/10000, round(IncAdv_Mat_Probit2(k+2)*10000)/10000, round(IncAdv_Mat_Probit4(k+2)*10000)/10000);
     fprintf(FID, '%s & %s & %s & %s \\\\ \n', blanks, strcat('(', num2str(IncAdv_Mat_Probit1_SE(k+2)), ')'), strcat('(', num2str(IncAdv_Mat_Probit2_SE(k+2)), ')'), strcat('(', num2str(IncAdv_Mat_Probit4_SE(k+2)), ')'));   
end

fprintf(FID, '\\hline \n');
             

fprintf(FID, '%s & %8.0f & %8.0f & %8.0f \\\\ \n ', 'Observations', length(VoteDem_1), length(VoteDem_1), length(VoteDem));
fprintf(FID, '%s & %8.3f & %8.3f & %8.3f \\\\ \n ', 'Log-likelihood', -fvalProbit1G, -fvalProbit2G, -fvalProbit4G);
fprintf(FID, '\\hline\\hline \n');
fprintf(FID, '\\multicolumn{3}{l}{\\tiny Bootstrap standard errors in parentheses are clustered at the district level.}\\\\ \n');
fprintf(FID, '\\end{tabular}\n');
fprintf(FID, '\\end{table} \n');
fclose(FID);

%%





Betas0 = vertcat(1, 1, zeros(length(RHS(1,:)), 1)); Theta_0 = vertcat(vertcat(0, 0), Betas0, 1);
[thetaTightProbit_InfoTest, fvalProbit_InfoTest] = fminunc(@(Theta) IndivVoterObjFunInfoTest_Senate(Theta,VoteDem, VoterIdeal, CandTypeDem, CandTypeRep, DemCandIdeal, DemPartyIdeal, RepCandIdeal, RepPartyIdeal, RHS, UtilFuncParam, ChoiceProbForm), Theta_0, optionsTight);

csvwrite('thetaTightProbit_InfoTest_Senate.csv', thetaTightProbit_InfoTest);
csvwrite('fvalProbit_InfoTest_Senate.csv', fvalProbit_InfoTest);


Betas0 = vertcat(1, 1, zeros(length(RHS(1,:)), 1)); Theta_0 = vertcat(vertcat(0, 0), Betas0, 0);
[thetaTightProbit_NoPolicy, fvalProbit_NoPolicy] = fminunc(@(Theta) IndivVoterObjFunNoPolicy_Senate(Theta,VoteDem, VoterIdeal, CandTypeDem, CandTypeRep, DemCandIdeal, DemPartyIdeal, RepCandIdeal, RepPartyIdeal, RHS, UtilFuncParam, ChoiceProbForm), Theta_0, optionsTight);

csvwrite('thetaTightProbit_NoPolicy_Senate.csv', thetaTightProbit_NoPolicy);
csvwrite('fvalProbit_NoPolicy_Senate.csv', fvalProbit_NoPolicy);


%%Creating Reciever Operator Chracteristic Graph

[PrProbit4] = PredictedProb_Senate(thetaTightProbit4G,VoteDem, VoterIdeal, CandTypeDem, CandTypeRep, DemCandIdeal, DemPartyIdeal, RepCandIdeal, RepPartyIdeal, RHS, UtilFuncParam, 1);
[FracIncorrectY0_Probit4, FracCorrectY1_Probit4]=get_ROC(PrProbit4, VoteDem);

[PrProbit_InfoTest] = PredictedProbInfoTest_Senate(thetaTightProbit_InfoTest,VoteDem, VoterIdeal, CandTypeDem, CandTypeRep, DemCandIdeal, DemPartyIdeal, RepCandIdeal, RepPartyIdeal, RHS, UtilFuncParam, 1);
[FracIncorrectY0_Probit_InfoTest, FracCorrectY1_Probit_InfoTest]=get_ROC(PrProbit_InfoTest, VoteDem);

[PrProbit_NoPolicy] = PredictedProbNoPolicy_Senate(thetaTightProbit_NoPolicy,VoteDem, VoterIdeal, CandTypeDem, CandTypeRep, DemCandIdeal, DemPartyIdeal, RepCandIdeal, RepPartyIdeal, RHS, UtilFuncParam, 1);
[FracIncorrectY0_NoPolicy, FracCorrectY1_NoPolicy]=get_ROC(PrProbit_NoPolicy, VoteDem);
%%
close 
myfiguresize =  [100 1 10 8];
set(gcf, 'PaperPosition', myfiguresize);
fh = figure(1);
get(gca)
set(gcf, 'color', 'white');

hold 'on'
line(FracCorrectY1_Probit4, FracIncorrectY0_Probit4)
line(FracCorrectY1_Probit_InfoTest, FracIncorrectY0_Probit_InfoTest, 'Marker', '+')
line(FracCorrectY1_NoPolicy, FracIncorrectY0_NoPolicy, 'Marker', 'o')
plot(0:0.05:1, 0:0.05:1, '--')
ylabel('True Positive Rate')
xlabel('False Positive Rate')
title('Reciever Operating Characteristics Curve')
legend('Full Model', 'Alternative Information Model', 'No Policy Model', 'No Predictive Power Model')
legend('Location', 'Southeast')
legend boxoff
hold 'off'


saveas(fh, 'ROC_Graphs_Senate', 'epsc'); 
close

matlabpool CLOSE